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Sk.Golam Ali, B.TalukdaiQ and S.K.Roy 

Department of Physics, Visva-Bharati University, Santiniketan 731235, India 

We study the dynamics of bright solitons in a Bose-Einstein condensate (BEC) confined in a 
highly asymmetric trap. While working within the framework of a variational approach we carry 
out the stability analysis of BEC solitons against collapse. When the number of atoms in the soliton 
exceeds a critical number N c , it undergoes the so called primary collapse. We find an analytical 
expression for N c in terms of appropriate experimental quantities that are used to produce and 
confine the condensate. We further demonstrate that, in the geometry of the problem considered, 
the width of the soliton varies inversely as the number of constituent atoms. 
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t^- ' 1. Introduction 

o : 

, A Bose-Einstein condensate (BEC) consists of trapped ultracold atoms all in the same quantum state. In this state 
the atoms lose their individual identities and behave as a single collective wave which is large enough to be optically 
imaged. In order to creat a BEC, atoms are first confined within a strong magnetic field and then the temperature of 
the atomic gas is continually lowered by laser and evaporative cooling until the condensate is formed. If one confines 
the BEC in only two directions, it will tend to disperse in the free direction. Because of the energetics involved, 
the atom-atom interaction in freely propagating BEC is characterized by the s-wave scattering length. The Feshbach 
resonance [lj allows one to continuously tune the scattering length from a positive to negative value (repulsive to 
attractive interaction ) by means of applied magnetic field. For attractive atomic interaction we can have coherently 
propagating matter-wave packets which travel over the BEC with neither attenuation nor change in shape. These are 
the so called bright solitons. For repulsive interaction we shall have dark solitons. Understandably, a bright soliton is 
a peak on the BEC while a dark soliton is a notch with a characteristic phase step across it. In a BEC of 7 Li atom 
i ^ , ' Rice and Paris teams Q produced bright solitons, each of which represents a condensate of actual atoms extracted 
from the main BEC. 

When the number of atoms in a bright soliton exceeds a critical value, it becomes unstable due to focusing nonlin- 
earity arising from the attractive atom-atom intarection. The transverse dimensions of the confinement then cause 
the soliton to collapse. This is often referred to as the primary collapse [s]. In the present work we shall envisage a 
variational study for the stability of bright soliton in a highly elongated trap and thereby calculate the critical number 
of atoms (N c ) that a soliton can hold before it undergoes the so called primary collapse. We shall see that the merit of 
our approach is its directness and simplicity because the variational method sought by us provides a straightforward 
\Q . analytical model to understand the dynamics of bright solitons. 

To extract the relevant physical information regarding stability and / or collapse we shall work within the framework 
of a mean field approximation. In this approximation the dynamics of a BEC is modelled by 3D Gross-Pitaevskii 
• i-H ' (GP) equation 
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where tp(f, t) is the macroscopic wave function of the condensate. This wave function is also called the order parameter. 
Here Uq — 47r ^ n as represents the interatomic interaction with a s , the two-particle s-wave scattering length and m, the 
mass of the atom. The wave function ip is normalized to the number of particles N in the condensates such that 



\ij)\ i dr = N. (2) 
The potential V(f) confines the atoms in a trap. For harmonic trapping, V(r) is given by 

V(r) = \mv 2 (\ 2 x x 2 + \ 2 y y 2 + \\z 2 ) . (3) 
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The parameters X Xl X y and A z describe anisotropy of the trap in the x, y and z direction respectively. We shall work 
with a highly asymmetric trap as determined by X x = A« = 1 and A z = - « 1. Here v z represents the frequency 
along the z directions and v = v r , the radial frequency. Our system of interest is thus a quasi-one dimensional (Q1D) 
BEC dispersing along the z direction. In the following we derive an appropriate version of the GP equation that will 
be useful to study the dynamics of BEC in highly asymmetric traps. 

We consider (1) in a geometry in which the trapping potential in z is much weaker than the corresponding potential 
in r = (x 2 + J/ 2 ) 5 - Further, we write the equation in terms of dimensionless variables defined by 



This gives 



z ^ 
r = ut, P=^-, s=—, ip(r,z,t) =u(p,s,T)/a$. (4) 



iu T + ^V 2 u - i (p 2 + \ 2 z s 2 ) u - ^^-\u\ 2 u = 0. (5) 
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Here ao = \ is the size of the ground state solution of the noninteracting GP equation. It is obvious that 

J\u\ 2 d 3 p = N. (6) 
We assume a separable ansatz for the solution of (5) such that Q 



u(p,s,t) = <P(p)£(s,t). (7) 



From (5) and (7) we have 



1 '-jV^+ipV), (8) 

where V 2 , stands for the Laplacian in the radial coordinate. In (8) the subscripts on £ stands for partial derivative with 

respect to that particular independent variable. More specifically, ^2s = fjf • This equation shows that the presence 
of atom-atom interaction does not permit clearcut separation of variables. However, the fourth term in equation (8) 
is quite small. Thus, 4> ma Y be assumed to satisfy 

+ \p 2 <t> = v p <l> (9) 

with v p being related to v r by a scale factor determined by the change of variables sought in (4). Equation (9) represents 
the well-known eigenvalue problem for the two dimensional harmonic oscillator with the ground sate solution given 
by 



Combining (8) and (9) we write 
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The low-frequency vibration along the z direction is quite unlikely to excite the two dimensional bosonic oscillator 
from its ground state. Thus (11) can be multiplied by </></>* and integrated over the p coordinate to get 



^ + \i2s-\\ls 2 t- 2 -^\i\ 2 i = Vp t (12) 
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Equation (12) represents the GP Equation for a Q1D trap. For a true ID system one does not expect collapse of the 
system with increasig number of atoms. But the use of Q1D trap in controlling the condensate motion may result in 
the collase of a BEC soliton when the number of atoms in it exceeds a critical value, say, N c . We shall work with 



(12) to provide an analytical model to study the collapse dynamics of bright solitons in a Q1D trap. Interestingly 
(12) can be written in a more convenient form by using the change of variable 



From (12) and (13) we get 



with 



as,r)= X (s,r)e- i ^. (13) 
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\x\ 2 ds = N/tt. (15) 
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Equation (14) represents the desired form of the evolution equation in which the atom-atom interaction is characterized 
by a negative scattering length. The realistic ID limit in (14) is not a true ID system because this equation involves 
the effect of transverse degrees of freedom through A z and ao. In this context we note that a similar equation with 
ip(r, z, t) chosen as u(p, s, r)/y/c^/N was used by Perez-Garcfa et al Q to qualitatively demonstrate that if the number 
of particles is large enough, the condensate is unstable and collapse occurs. We are , however, interested to derive a 
straightforward analytical model to understand the collapse dynamics and thereby provide a quantitative estimate 
for N c . To that end we convert, in section 2, the initial-boundary value problem in (14) to a variational problem. 
In particular, we present an expression for the Lagrangian density and a trial wave function involving variational 
parameters to study the stability of bright solitons against collapse. We also obtain the evolution equations for these 
parameters. In section 3 we judiciously use the derived evolution laws to study the soliton dynamics with particular 
emphasis on the stability of solitons against collapse. 



2. Variational formulation 

The action principle 
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with the Lagrangian density given by 



c = | (xx* - x*x T ) + \>? z s\ x * + — xV 2 + \xlXs, (17) 
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is equivalent to (14). We shall use this expression for C to studythe dynamics of bright solitons in terms of a 
variational method often called the Ritz optimization procedure [5J. In this procedure the first variation of the 
variational functional is made to vanish within a set of suitable chosen trial function such that the field theoretical 
problem under consideration reduces to a simple problem of point mechanics. For the negative scattering length the 
Gaussian trial function for x(s, r) is a very reasonable ansatz. Thus we write 

x (s, t) = A(t) exp [-s 2 /2a 2 (r) + i6(r)s 2 /2] ■ (18) 

Here A(t) is a complex amplitude, a(r) - the width of the distribution and 6(r) - the frequency chirp. The phase of 
the condensate S(t) is defined by A(t) — \ A(T)\e s ^ T \ The amplitude A(t), width a(r) and the chirp b(r) will all vary 
with the time parameter r. The initial condensate at rest will have da(T)/dr = 0. 

Inserting the trial function in (18) into the variational principle stated in (16) we obtain a reduced variational 
problem 

5 J (C)dr = (19) 

with 
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Here Cq represents the result obtained by inserting the Gaussian ansatz (18) into the Lagrangian density in (17). It 
is rather straight forward to perform the integration in (20) and get 

(£) = ypH 1 - (AA* - A* A T ) a + h T a 3 AA* + 
\\ 2 a 3 AA* + -^=—aA 2 A* 2 + U -+b 2 a 3 ) AA*) . (21) 



Let us now obtain the variational equations for the Gaussian parameters A(t), A*(t), a(r) and b(r) which follow 
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from the vanishing conditions of ^nr, xzr-, and . These equations are given by 
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and 
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From (22) and (23) we have found 



5 {C) - ha 3 AA* -~ (a 3 AA*) = 0. (25) 



such that 



4- (aAA*) = (26) 
dr 



a\A\ 2 — Q, a constant. (27) 



The constant Q is simply related to the number of particles in the condensate since the value of the integral (15) is 
y/wa\A\ 2 . Combining (25) and (27) we get 

b=4~ (In a). (28) 
dr 

Equations (27) and (28) clearly show that if we can derive a method to calculate the values of a(r), the other parameters 
of the condensate will be automatically determined. Fortunately, (22), (23), (24) and (28) can be combined to write 
a second-order ordinary differential equation, the first integral of which gives 

1 / da\ 2 1 , 2 2 [2Na s 1 1 
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with E, the constant of integration. 



The equation for a(r) in (29) is related to the motion of a particle in a potential field V(a) so that 



f /da N 



2{Tr) + V (a) = E. (30) 
Here 

V(a) = \\W + - + A ' P = \^~- ( 31 ) 
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Thus one would like to interpret the constant of the motion E as the total energy of the particle. It is easy to 
verify that there is no physical uncertainty in the identification sought because H-Gds represents the left side 
of (29). Here Tic stands for the result of the Hamiltonian density calculated from (18) and rewritten by using the 
Gaussian ansatz. Obviously E is determined by the initial conditions of the second order differential equation from 
which (30) has been extracted. It is not difficult to solve (30) and look for the dynamics of the condensate. However, 
the analysis of the equilibrium point obtained from the extrimum of V(a) written as 

— r- 2 = 32) 

da 

can give some illuminating results. 



3. Dynamics of bright solitons 

For bright solitons the nonlinear interaction is attractive and the scattering length a s < 0. In this case we shall use 
P = —\P\ and carry out the subsequent analysis by using only the numerical values of a s . We shall make use of (32) 
to derive a simple physical picture for the collapse dynamics of bright solitons when the trap of the BEC is relaxed 
in one direction. From (31) and (32) with P = — \P\ we get 

X 2 z a 4 + \P\a- 1 = 0. (33) 

The equilibrium point determined by (33) should be a minimum for (14) to support a soliton solution. This gives 

\ 2 z a 4 -2\P\a + 3 = f3a 4 , (3 >0. (34) 

Eliminating \P\ from (33) and (34) we find that 

a =G^) 4 '" >3A - (35) 

is a particular solution of (33) and (34). From (35) and (33) or (34) we get 



(36) 



The form of (36) imposes a further restriction on the values of (3 than that given in (35) and sets a lower bound for 
it. Using (3 — "f\ 2 z we write (36) in the form 

W = ^4^ (37) 
(7-3)3 

Thus non zero values of P will be obtained for 7 > 4 only. For 7 = 4 the interaction term vanishes and GP equation 
becomes linear and soliton formation becomes impossible. 
From (31) and (37) we obtain an expression 

V 2 \a s \ ( 7 _ 3)4 



for the number of atoms in the quasi-lD soliton. 
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FIG. 1: The potential V(a) as a function of a for \ z = 

In Fig. 1 we plot the potential V(a) in (31) as a function of a for a s — —1.59 x 10 _4 /im, the scattering length of 7 Li 
as used in the experiment of Strecker et al. 0]. We have chosen to work with A z = jgg. In this figure we have four 
curves represented by ^4(0), V§{a), Vq{o) and Vj(a) corresponding to 7 = 4, ,5, 6, and 7 repectively. A common 
feature of all these potentials is that each of them exhibits a minimum. The curve for Vj(a) represents a potential 
well between a\ — 5.0114 and a<i — 10.5468. The minimum of the well is negative. A mechanical analogy suggests 
a solution which oscillates between the zeros of V7 (a) . In this case, spreading of the BEC is stopped at a — ai by 
nonlinear effects which subsequently compress the BEC back to the initial width. This behavior is repeated in an 
oscillatory manner. In this situation the BEC soliton will become unstable and lead to a mechanical collapse [7]. A 
similar situation arises for other values of 7 > 6. 

For 7 = 6 the potential well degerates into a single point such that Ve(a) touches the a axis at a particular point 
where the potential has a stable minimum. Understandably, a particle releasesd at this point will stay there. In the 
present context this implies that for our chosen value of A z and 7 = 6 the BEC bright soliton will be critically stable. 
Using 7 = 6 we get 

N c = 0.8774, /J^VaT- (39) 
V 2 \a„\ 

After the number of atoms exceeds this critical number the soliton becomes unstable. More than a decade ago 
Ruprecht et al. [8[ used a purely numerical routine to set a limit for the critical number of atoms after which the BEC 
with attractive two-body interaction becomes unstable. This limit was examined by Gammal et al. Q for different 
trapping geometries. Interestingly, the analytical expression in (39) is in agreement with the observations of refs. 

8 and 9. From (39) it is clear that ^"J" 3 ' <§; 1. This represents the well-known relation for the existance of stable 
solitons fioj |. 

It will be interesting to see what happens if 7 < 6. In order that we look at the curve represented by Vi(a). For 
7 = 4, \P\ = and we do not have a nonlinear term in the GP equation. In this case, no soliton can be formed. The 
minimum of V4 is positive. This appears to suggest that weaker nonlinearies leading to potential curves with positive 
minima will not be able to produce matter- wave bright solitons. In our figure Vs(a) represents one such curve. 

Equations (35) and (37) can be combined to write 



o(r) = ^, (40) 

where 
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(41) 



The result in (40) is remarkable and implies that, at a given instant of time, if the number of atoms in the soliton 
increases, it becomes narrower. It is a real curiosity to note that the result in (40) holds good even in the absence of 
trapping [H 
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